home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Libris Britannia 4
/
science library(b).zip
/
science library(b)
/
SCIENTIF
/
0807.ZIP
/
GALILEAN.BAS
< prev
next >
Wrap
BASIC Source File
|
1987-09-30
|
13KB
|
463 lines
' Program "GALILEAN"
' copyright (C) 1986 by David Eagle
' public domain for IBM-PC on November 23, 1986
' IBM-PC << QuickBASIC Compiler Version 3.0 >>
' determines position of Jupiter's great satellites
' x-position (positive west of Jupiter; units of Jupiter radii)
' y-position (positive north of Jupiter; units of Jupiter radii)
' position angle (relative to inferior junction with Jupiter; degrees)
' Jupiter semi-diameter (units of arc seconds)
' Julian Date
'********************************************************************
OPTION BASE 1
DEFDBL A-Z
DIM SHARED MONTH$(12),SATELLITE$(4),R(4),X(4),Y(4),U(4)
COMMON SHARED JULIAN.DATE0,HOUR%,MINUTE%,LCT$,MONTH%,DAY%,YEAR%
COMMON SHARED DST%,ZONE%,JULIAN.DATE,SEMIDIAMETER,CDATE$,TDELAY
CONST PI=3.141592653589793D0
CONST PIDIV2=.5D0*PI
CONST PI2=2D0*PI
CONST RTD=180D0/PI
DEF FNMOD(X)=X-PI2*INT(X/PI2)
DEF FNJDATE(MONTH%,DAY%,YEAR%)
'JULIAN Date function
F=367D0*YEAR%-INT(7*((YEAR%+INT((MONTH%+9)/12))/4))
F=F+INT(275*MONTH%/9)+DAY%+1721028.5D0
FNJDATE=F-INT(3*(INT((YEAR%+SGN(MONTH%-9)*INT(ABS((MONTH%-9)/7)))/100)+1)/4)
END DEF
DEF FNASN(X)
' inverse sine function
IF ABS(X)>=1D0 THEN
FNASN=SGN(SINANGLE)*PIDIV2
ELSE
FNASN=ATN(X/SQR(1D0-X^2))
END IF
END DEF
' calendar months
MONTH$(1)="January"
MONTH$(2)="February"
MONTH$(3)="March"
MONTH$(4)="April"
MONTH$(5)="May"
MONTH$(6)="June"
MONTH$(7)="July"
MONTH$(8)="August"
MONTH$(9)="September"
MONTH$(10)="October"
MONTH$(11)="November"
MONTH$(12)="December"
' satellite names
SATELLITE$(1)="Io"
SATELLITE$(2)="Europa"
SATELLITE$(3)="Ganymede"
SATELLITE$(4)="Callisto"
TDELAY=10
'********************************************************************
CLS
PRINT
PRINT "Program Galilean"
PRINT "(C) Copyright 1986 by David Eagle"
PRINT
PRINT "Microsoft QuickBasic Compiler"
PRINT "(C) Copyright Microsoft Corp. 1982-1987"
PRINT
CALL KEYCHECK
CLS
PRINT
PRINT
INPUT "Program introduction ( y = yes, n = no ) ";INTRO$
IF INSTR("yY",INTRO$) THEN CALL INTRODUCTION
DO
CLS
PRINT
PRINT
PRINT
PRINT "Calendar date ( month [ 1 - 12 ], day [ 1 - 31 ], year [ YYYY ] )"
PRINT "< NOTE: B.C. dates are negative, A.D. dates are positive >"
PRINT "< For example, October 21, 1986 is input as 10,21,1986 >"
INPUT MONTH%,DAY%,YEAR%
PRINT
PRINT "Local civil time ( hours [ 0 - 23 ], minutes [ 0 - 59 ] )"
PRINT "< For example, 8:30 pm is input as 20,30 >"
INPUT HOUR%,MINUTE%
PRINT
PRINT "Time zone ( 0 - 23 )"
PRINT "< For example, Mountain Standard Time (MST) is time zone 7 >"
INPUT ZONE%
PRINT
PRINT "Daylight Savings Time ( y = yes, n = no )"
INPUT DST.FLAG$
IF INSTR("yY",DST.FLAG$) THEN
DST%=1
ELSE
DST%=0
END IF
JULIAN.DATE0=FNJDATE(MONTH%,DAY%,YEAR%)
CALL MAIN.DRIVER
CLS
PRINT
PRINT
PRINT "Another selection ( y = yes, n = no )"
INPUT SELECTION$
LOOP UNTIL INSTR("nN",SELECTION$)
END
'********************************************************************
SUB MAIN.DRIVER STATIC
' main driver subroutine
FIRSTPASS$="TRUE"
RESPONSE%=0
JULIAN.DATE0=JULIAN.DATE0-1D0
LCT$=STR$(HOUR%)+" hours "+STR$(MINUTE%)+" minutes"
WHILE RESPONSE%<>4
RESPONSE%=0
WHILE RESPONSE%<=0 OR RESPONSE%>5
CLS
CALL SELECTION.MENU(RESPONSE%)
WEND
IF RESPONSE%=3 OR FIRSTPASS$="TRUE" THEN
JULIAN.DATE0=JULIAN.DATE0+1D0
CALL GREGORIAN.DATE(JULIAN.DATE0+.5D0)
CALL SAT.POSITIONS
END IF
IF RESPONSE%=1 THEN CALL DISPLAY.DATA
IF RESPONSE%=2 THEN CALL GRAPHICS
FIRSTPASS$="FALSE"
WEND
END SUB
'********************************************************************
SUB ELAPSED.DAYS(EDAYS) STATIC
' elapsed days subroutine
GMT=HOUR%+MINUTE%/60D0-DST%+ZONE%
A=(JULIAN.DATE0-2415020D0)/36525D0
JULIAN.DATE=JULIAN.DATE0+GMT/24D0+(.41D0+1.2053D0*A+.4992D0*A^2)/1440D0
EDAYS=JULIAN.DATE-2415020D0
END SUB
'********************************************************************
SUB SAT.POSITIONS STATIC
' satellite positions subroutine
CALL ELAPSED.DAYS(D2)
V=FNMOD(2.34974D0+.0000194756D0*D2)
M=FNMOD(6.256586D0+.01720197D0*D2)
N=FNMOD(3.93272126D0+.0014501127D0*D2+.00576D0*SIN(V))
J=FNMOD(3.8684699D0+.015751909D0*D2-.00576D0*SIN(V))
O=FNMOD(.0334405D0*SIN(M)+.000349D0*SIN(2D0*M))
P=FNMOD(.0969007D0*SIN(N)+.0029147D0*SIN(2D0*N))
K=FNMOD(J+O-P)
R1=1.00014D0-.01672D0*COS(M)-.0014D0*COS(2D0*M)
R2=5.20867D0-.25192D0*COS(N)-.0061D0*COS(2D0*N)
D3=SQR(R1^2+R2^2-2D0*R1*R2*COS(K))
S4=FNASN(.00047726615D0/D3)
SEMIDIAMETER=206264.806D0*S4
S2=R1*SIN(K)/D3
S3=FNASN(S2)
L1=4.154756D0+.0014502D0*D2+.00576D0*SIN(V)+P
A=.05358D0*SIN(L1+.77667D0)
D4=A-.03752D0*S2*COS(L1+.4189D0)-.02286D0*((R2-D3)/D3)*SIN(L1-1.73486D0)
A=D2-D3/173.1446D0
U(1)=FNMOD(1.475686D0+3.550102027D0*A+S3-P)
U(2)=FNMOD(.724338D0+1.767872488D0*A+S3-P)
U(3)=FNMOD(1.9194607D0+.876757718D0*A+S3-P)
U(4)=FNMOD(3.07803823D0+.375036004D0*A+S3-P)
G=FNMOD(3.269D0+.87808691D0*A)
H=FNMOD(5.4297D0+.376454063D0*A)
R(1)=5.9061D0-.0244D0*COS(2D0*(U(1)-U(2)))
R(2)=9.3972D0-.0889D0*COS(2D0*(U(2)-U(3)))
R(3)=14.9894D0-.0227D0*COS(G)
R(4)=26.3649D0-.1944D0*COS(H)
U(1)=FNMOD(U(1)+.0082296D0*SIN(2D0*(U(1)-U(2))))
U(2)=FNMOD(U(2)+.018727D0*SIN(2D0*(U(2)-U(3))))
U(3)=FNMOD(U(3)+.003037D0*SIN(G))
U(4)=FNMOD(U(4)+.014748D0*SIN(H))
FOR I%=1 TO 4
X(I%)=R(I%)*SIN(U(I%))
Y(I%)=-R(I%)*COS(U(I%))*SIN(D4)
NEXT I%
END SUB
'********************************************************************
SUB GREGORIAN.DATE(JDATE) STATIC
' Gregorian Date subroutine
IF JDATE<2299161D0 THEN
A=JDATE
ELSE
A=INT((JDATE-1867216.25D0)/36524.25)
A=JDATE+A-INT(A/4)+1D0
END IF
B=A+1524D0
C=INT((B-122.1)/365.25D0)
D=INT(365.25D0*C)
E=INT((B-D)/30.6001D0)
DAY%=B-D-INT(30.6001D0*E)
IF E<13.5 THEN
MONTH%=E-1
ELSE
MONTH%=E-13
END IF
IF MONTH%>2.5 THEN
YEAR%=C-4716
ELSE
YEAR%=C-4715
END IF
CDATE$=MONTH$(MONTH%)+STR$(DAY%)+","+STR$(YEAR%)
END SUB
'********************************************************************
SUB SELECTION.MENU(SELECTION%) STATIC
' selection menu subroutine
CLS
PRINT
PRINT
PRINT " *** Please Note ***"
PRINT
PRINT
PRINT " After a data or graphics display page appears on the screen, it may"
PRINT " be saved to a printer by pressing the <Shift> <PrtSc> key combination."
PRINT
PRINT " Printing a graphics screen requires that the program 'GRAPHICS.COM'"
PRINT " from the system disk be executed before running this program."
PRINT
PRINT " The program will automatically display data or graphics screens every"
PRINT TDELAY;"seconds. The user may manually cycle the screen by pressing any key."
PRINT
PRINT
PRINT
PRINT "Would you like to change the screen display interval ( y = yes, n = no )"
INPUT A$
IF INSTR("yY",A$) THEN
PRINT
PRINT
INPUT "New time interval ( seconds )";TDELAY
END IF
CLS
PRINT
PRINT
PRINT TAB(25);"Galilean Menu"
PRINT
PRINT
PRINT TAB(20);"< 1 > Display Data"
PRINT
PRINT TAB(20);"< 2 > Graphics"
PRINT
PRINT TAB(20);"< 3 > Continue"
PRINT
PRINT TAB(20);"< 4 > End"
PRINT
PRINT
PRINT "Selection"
INPUT SELECTION%
END SUB
'********************************************************************
SUB DISPLAY.DATA STATIC
' display data subroutine
CLS
PRINT
PRINT TAB(5);"Calendar date";TAB(60-LEN(CDATE$));CDATE$
PRINT
PRINT TAB(5);"Local civil time";TAB(60-LEN(LCT$));LCT$
PRINT
PRINT TAB(5);"Julian Date";
PRINT TAB(48);
PRINT USING "########.###";JULIAN.DATE
PRINT
PRINT TAB(5);"Semidiameter ( arc seconds )";
PRINT TAB(48);
PRINT USING "########.###";SEMIDIAMETER
PRINT
PRINT
PRINT TAB(5);"Satellite X-position Y-position Angle"
PRINT TAB(5);"--------- ---------- ---------- -----"
FOR I%=1 TO 4
A$=STR$(INT(U(I%)*RTD))
PRINT
PRINT TAB(5);SATELLITE$(I%);
PRINT TAB(22);
PRINT USING "###.###";X(I%);
PRINT TAB(39);
PRINT USING "###.###";Y(I%);
PRINT TAB(59-LEN(A$));A$
NEXT I%
TDELAY=20
CALL KEYCHECK
END SUB
'********************************************************************
SUB GRAPHICS STATIC
' graphics subroutine
SCREEN 2
PRINT
PRINT TAB(5);CDATE$;
PRINT TAB(75-LEN(LCT$));LCT$
PRINT TAB(38);"North"
LOCATE 9,1
PRINT TAB(70);"West"
A$="JUPITER"
LOCATE 5,1
FOR J%=1 TO 7
PRINT TAB(40);MID$(A$,J%,1)
NEXT J%
CIRCLE (315,110),10
PAINT (315,110)
LOCATE 16,1
FOR I%=1 TO 4
X=315+10*X(I%)
Y=110-4.4444*Y(I%)
S$="^ "+SATELLITE$(I%)
CIRCLE (X,Y),5
PAINT (X,Y)
IF (X+LEN(S$))>550 THEN
S$=SATELLITE$(I%)+" ^"
X=X-7*LEN(S$)
END IF
PRINT TAB(X/8);S$
NEXT I%
A=TIMER
B=TIMER+TDELAY
A$=""
WHILE A$="" AND A<B
A=TIMER
A$=INKEY$
WEND
SCREEN 0
END SUB
'********************************************************************
SUB INTRODUCTION STATIC
' program introduction subroutine
CLS
TDELAY=30
PRINT
PRINT
PRINT TAB(4);"GALILEAN is an interactive QuickBASIC program which can be"
PRINT TAB(4);"used to determine the position of the Galilean satellites"
PRINT TAB(4);"relative to Jupiter. This information is useful for such"
PRINT TAB(4);"activities as astronomical observations, astrophotography"
PRINT TAB(4);"and the study of occultations between the moons."
PRINT
PRINT TAB(4);"This program provides the position of each of the four"
PRINT TAB(4);"Galilean moons in the units of the radius of Jupiter. The"
PRINT TAB(4);"x-position of each satellite is measured positive west of"
PRINT TAB(4);"Jupiter and the y-position is measured positive north."
PRINT TAB(4);"It also calculates the 'position angle' of each satellite."
PRINT TAB(4);"This angle helps determine if any of the satellites are"
PRINT TAB(4);"in front of or behind Jupiter. A position angle near 0"
PRINT TAB(4);"or 360 degrees is called 'inferior conjunction' and an"
PRINT TAB(4);"angle near 180 degrees is called 'superior conjunction'."
PRINT TAB(4);"A satellite is in front of Jupiter at inferior conjunction"
PRINT TAB(4);"and behind Jupiter at superior conjunction."
CALL KEYCHECK
CLS
PRINT
PRINT
PRINT TAB(4);"Inputs required by program 'GALILEAN' include the user's"
PRINT TAB(4);"observation date, the local time and time zone. The user"
PRINT TAB(4);"must also specify if Daylight Savings Time is in effect."
PRINT TAB(4);"The date is input in numerical format and the local time"
PRINT TAB(4);"in 24 hour format. The time zone will be an integer number"
PRINT TAB(4);"between 0 and 23."
PRINT
PRINT TAB(4);"The 'Galilean Menu' will allow the user to display the"
PRINT TAB(4);"data and/or graphics for the selected date. The 'Continue'"
PRINT TAB(4);"selection of this menu will calculate the information for"
PRINT TAB(4);"the next day. The user can then choose to display data or"
PRINT TAB(4);"and/or graphics for the new date. Please note that you can"
PRINT TAB(4);"print the data or graphics screen by pressing the <Shift>"
PRINT TAB(4);"<PrtSc> keys. The 'End' selection of the 'Galilean Menu'"
PRINT TAB(4);"will allow you to end the current session and select a new"
PRINT TAB(4);"date and/or time."
CALL KEYCHECK
END SUB
'********************************************************************
SUB KEYCHECK STATIC
' check user response subroutine
PRINT
PRINT
PRINT TAB(20);"< press any key to continue >"
A=TIMER
B=TIMER+TDELAY
A$=""
WHILE A$="" AND A<B
A=TIMER
A$=INKEY$
WEND
END SUB